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ABSTRACT 


A  simplified  but  nonlinear  model  is  used  for  the  autopilot  in  the 
single-axis  missile  interceptor  problem.  A  recent  exact  design  technique 
for  uncertain  nonlinear  feedback  systems,  is  used  to  design  the  autopilot. 

A  large  number  of  target  initial  range,  position,  velocity  and  acceleration 
values  were  used.  Also,  uncertainty  in  the  aerodynamic  missile  parameters 
was  allowed.  The  objective  was  to  guarantee  that  the  autopilot  response 
satisfied  specified  tolerances  over  the  above  ranges.  This  was 
satisfactorily  achieved.  This  autopilot  design  was  then  used  in  the 
simulation  of  the  single-axis  interceptor  system,  over  the  above  set  of 
target  trajectories.  The  results  were  highly  satisfactory. 
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APPLICATION  OF  SYNTHESIS  TECHNIQUE  FOR  NONLINEAR  UNCERTAIN 
SYSTEMS  TO  FLIGHT  CONTROL  OF  ADVANCED  MISSILES 

1.  INTRODUCTION 

1.1.  In  many  control  systems  there  is  a  given  constrained  part,  denoted 
here  as  the  Plant,  whose  output  is  the  system  output  of  interest.  Math¬ 
ematically,  the  plant  is  represented  by  nonlinear  differential  equations. 
Often,  the  parameters  of  the  differential  equations  are  not  precisely  known 
and/or  may  change  in  time.  Also,  there  may  be  external  disturbances,  such  as 
wind  gusts,  not  known  in  advance.  Generally,  it  is  desired  that  the  system 
output  satisfy  given  performance  specifications,  despite  the  parameter  un¬ 
certainty  and  the  external  disturbances,  so  a  feedback  structure  is  necessary 

There  have  not  been  available  exact  techniques  for  designing  feedback 
systems  around  nonlinear  uncertain  plants,  so  a  variety  of  approximate 
techniques  have  emerged.  These  are  usually  supplemented  by  considerable 
testing  and  experimentation.  However,  there  has  recently  appeared  an  exact 
synthesis  procedure  for  a  large  class  of  nonlinear  plant  (even  nonlinear  time 
varying)  with  significant  parameter  uncertainty  [1,2].  The  specialty  of  this 
technique  Is  synthesis  to  satisfy  assigned  performance  specifications.  In 
fact,  the  assignment  of  such  performance  specifications  is  essential  —  which 
is  generally  not  an  unreasonable  demand.  But  even  in  those  cases  where  there 
is  uncertainty  as  to  what  constitutes  reasonable  performance,  one  may  use  his 
best  judgement  to  formulate  them  anyhow  and  proceed  with  the  design  technique 
If  the  performance  specifications  are  unreasonable,  the  design  which  emerges 
will  reveal  this,  typically  by  signal  levels  at  some  points  which  exceed  the 
saturation  values  of  plant  elements.  The  specifications  can  then  be  suitably 
modified. 

1.2.  The  Nonlinear  Design  Philosophy 

To  simplify  the  presentation,  consider  a  nonlinear,  t ime- invariant  plant 
set  W  ■  { w)  (a  set  because  of  uncertainty).  Each  parameter  combination 
gives  a  different  w  .  The  basic  idea  is  to  convert  this  set  into  a  linear 
time- Invariant  (denoted  by  fcti  )  plant  set  P  =»  { p(s) )  ,  such  that  P  _i_s. 
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precisely  equivalent  to  the  nonlinear  set  W  providing  the  output  of  any 
P(s)  In  P  belongs  to  the  set  of  acceptable  outputs.  To  help  grasp  this 
idea,  consider  a  black  box  Bj  which  contains  a  nonlinear  element  Wq  . 

Let  the  input  be  Xj  and  the  resulting  output  be  y^  -  wQ(xj)  .  Now  one 
could  (in  most  cases)  build  another  black  box  Bg  containing  a  linear 
element  p  such  that  for  this  same  input  x^  the  output  is  y1  -  p(x^)  . 

If  you  are  not  allowed  to  apply  any  other  input  except  x^  ,  it  is  impossible 
for  you  to  determine  which  is  the  linear  and  which  the  nonlinear  box.  The 
following  is  a  simple  technique  for  finding  p  .  Obtain  the  Laplace  trans¬ 
forms  of  y,  (t)  ,  x j  ( t )  :  Yj  (s)  -  l-y f  ( t )  ,  X,(s)  -  Lx^t)  and  let  P(s)  - 
Y i  (s) /X |  (s)  ■  Then  p(t)  -  L  1 P (s)  and  y^t)  -  p(t)*Xj(t)  ,  where  * 
indicates  linear  convolution  of  signals.  Hence  two  necessary  conditions  for 
this  equivalence  of  iti  p  to  nonlinear  w^  are: 

1.  The  output  yj  ■  Wq(xj)  is  unique  or  equivalently  given  y^  ,  the  input 
x,  ■  wQ' (yj)  is  unique. 

2.  The  signals  y^t)  ,  Xj(t)  are  Laplace  transformable.  This  condition 
is  difficult  to  violate. 

Condition  1  excludes  a  large  nonlinear  class,  e.g.  hard  saturation. 
Nevertheless,  this  technique  has  been  applied  to  the  latter  case,  by  replacing 
the  hard  saturation  by  very  low  gain  over  the  applicable  intervals  [2].  But 
obviously,  it  is  then  essential  that  the  specified  set  of  desired  system  out¬ 
puts,  be  consistent  with  the  limitation  due  to  saturation.  Actually,  a  good 
argument  can  be  made  for  the  violation  of  Condition  1,  providing  (a)  the 
single  output  y^  can  be  associated  with  a  compact  set  X^  =  {x}  ,  i.e.  a 
set  which  can  be  approximated  as  accurately  as  desired  by  a  finite  number  of 
elements  and  (b)  a  small  enough  change  in  output  y^  results  in  a  small 
change  in  the  associated  x^  .  But  this  has  not  as  yet  been  rigorously  proven. 

In  any  case,  in  the  above,  only  one  nonlinear  plant  Wq  and  one  output 
y I ( t )  were  considered.  The  idea  can  be  extended  to  uncountable  sets 
W  *  {w}  and  V  *  (y(t))  ,  as  follows.  We  work  backwards  from  the  output  of 
the  nonlinear  plant,  because  in  a  synthesis  problem  one  generally  knows  what 
kind  of  outputs  he  would  like  to  have  in  response  to  the  command  input  (not 
the  plant  input  —  for  example  the  latter  could  be  the  elevator  control  surface 
while  the  former  is  the  pilot  stick  motion).  Due  to  uncertainty,  there  is  a 
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s«t  of  nonlinear  plants  W  •  {w}  .  Take  any  pair  y*  €  V  ,  Wj  €  III  and 

find  xi(t)  ■  w,'(y*(t))  .  Then  find  pi  ,  the  ttl  equivalent  of  w, 

J  J  1  J  |  J  J 

(with  respect  to  y  )•  es  shown  above  e.g»9  by  finding  Y  (s)  *  Ly  (t)  , 

Xj (s)  ■  Lxj ( t )  and  setting  Y*(s)/xj(s)  ■  fj  (s)  ■  Lpj(t)  •  This  Is 

repeated  over  the  Wj  €  W  and  the  y*  €  V  giving  a  set  P  ■  {fj)  .  The 

set  P  J_s  the  ttl  equivalent  of  W  ,  with  respect  to  the  set  V_  .  1  f  W 

has  10  elements  and  V  has  20,  then  ?  has  200  elements.  In  practice,  both 

V  and  y  are  In  general  Infinite  and  so  is  P  .  Of  course,  in  the  execution 
of  the  design  technique  a  finite  number  of  V  and  Y  elements  are  taken. 

The  above  tti  equivalence  idea  is  the  heart  of  the  method.  Once  we 
have  a  linear  time- invariant  plant  set  P  ,  we  have  a  purely  tti  problem. 

It  is  essential  that  we  be  able  to  solve  this  tti  problem.  For  only  then 
it  is  guaranteed  that  the  actual  output  is  indeed  a  member  of  the  acceptable 
set  y  ,  no  matter  which  w  €  W  happens  to  be  the  actual  plant.  Thus,  the 
solution  (  F(s)  ,  G(s)  in  Figure  1  if  tti  compensation  is  used),  to  the 
tti  problem  is  guaranteed  to  also  be  the  solution  of  the  nonlinear  un¬ 
certainty  problem.  It  is  important  to  emphasize  that  there  are  no  approx¬ 
imations  involved  here.  The  procedure  is  exact  and  theoretically  rigorous. 

In  the  above,  only  a  single  command  input  r-r  in  Figure  1  was 

d 

assumed,  and  the  output  set  Y  is  the  set  of  acceptable  responses  to  r  (t) 
i.e.  in  Figure  I  it  is  required  that  c(t)  €  V  .  This  is  a  highly  exceptional 
case.  Usually  one  has  a  set  R  ■  { r )  of  typical  inputs.  Suppose  one  wants 
the  closed-loop  system  to  essentially  behave  like  a  tti  one,  in  response  to 
these  inputs,  with  transfer  function  T(s)  f  T  ,  the  set  of  acceptable  transfer 
functions.  Note  that  a  set  of  acceptable  l T ( s) >  must  be  specified  when  there 
is  plant  uncertainty,  because  invariance  is  impossible  to  achieve.  Then  the 
overall  total  set  V  •  (Y  I  is  obtained  by  generating  Yb(s)  ■  R  (s)T^(s) 
for  the  R  €  R  ,  i  €  T  .  One  then  finds  the  tti  equivalent  of  this  total 

V  set. 

It  is  worth  emphasizing  that  the  synthesis  technique  has  even  greater 
flexibility.  Suppose  an  overall  closed-loop  nonlinear  system  is  desired,  in 
the  following  sense:  For  input  r^  one  wants  system  response  characterized 
by  •  (Ta(s) I  ,  but  for  input  r^  one  wants  system  response  characterized 
by  set  T^  ■  (T^(s)l  .  For  example,  suppose  r^  is  a  unit  step  and  its 


FIGURE  1  Nonlinear  set  W  is  replaced  by  lti  equivalent  set 


FIGURE  3  Missile-target  geometry 


acceptable  response  function  set  is  T  .  Suppose  r^  is  a  unit  ramp  but 
ona  does  not  want  the  system  response  to  be  the  integral  of  its  response 
when  the  input  was  a  unit  step,  but  of  a  different  nature  entirely.  This 
can  be  done.  The  F  in  Figure  1  must  then  be  nonlinear.  See  Reference  2 
for  an  example  of  such  design.  The  output  c(t)  due  to  a  two-value  step 
need  not  be  twice  that  of  a  unit-valued  step  etc.,  if  one  does  not  want  it  so. 

It  Is  worth  reemphasizing  the  importance  of  the  synthesis  approach.  One 
must  formulate  the  set  of  desired  plant  outputs  or  alternatively  its  desired 
set  of  inputs.  Otherwise  it  is  not  possible  to  find  the  tt)  equivalent  set 
P  ,  which  is  the  basis  of  the  technique.  It  is  also  worth  Including  here,  as 
Appendix  I,  the  design  results  for  an  application  of  this  design  technique, 
recently  done  for  the  U.S.  Air  Force.  Here,  the  set  of  acceptable  plant  out¬ 
puts  was  very  clearly  apriori  specified  at  the  very  beginning.  Unfortunately, 
in  our  present  case,  however,  a  great  deal  of  work  was  required  to  obtain  the 
set  of  acceptable  outputs. 


2 .  STATEMENT  OF  THE  PROBLEM 

2.1.  The  problem  considered  is  the  single-axis  missile  interceptor 
problem,  shown  in  Figure  2.  Feedback  around  the  body  pitch-rate  is  usually 
used.  The  nonlinear  missile  equations  of  motion  are  usually  linearized. 

But  If  the  missile  has  a  large  flight  envelope,  then  such  linearization  may 
be  inadequate.  The  objective  here  is  to  use  nonlinear  equations  of  motion 
and  the  ttl  plant-equivalent  method  described  in  Section  1,  to  perform  the 
design  of  the  auto-pllct  in  Figure  2. 

The  equations  describing  Figure  2,  are 
Autopl lot  (Figure  3). 

a  •  angle  of  attack  (radians),  9  ■  q  ■  body  axis  pitch  rate  (r/s,  s» seconds), 
\  •  line  of  sight  angle  (r),  y  ■  velocity  angle  tangent  to  trajectory  (r) , 

6  •  elevator  deflection  (r). 
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Y  -  C^I^Wa  ♦  ^(4)4) 
e  -  CjlK^cOa  ♦  *^(4)4  ♦  C2$] 
o  •  0  -  y 


(U-c) 


V 

9  " 
N  - 

\a 


2 

1  ®^M^m  Per  **c‘  •  ^3  "  QD/A  per  sec  .  , 

9.81  m/s^  (meters/sec*)  ,  Q  ■  176  kgm 

70  k  (kgm)  ,  V  -  200  m/s  ,  0  -  .17 

w 

-  CU(17+ 1.3a)  ,  *L4  ■  CL6(15+ 1.1«) 

-  W-»  ‘  -  *»>  • 


} 


K^Dm  sec. 

K  -  -1000 

mq 

A  -  3.8 


(2) 


(3a-d) 


In  (3a-d),  the  coefficients  are  Independently  uncertain. 


CL«'CU-C™,’C«  €  '0.8.  1.21  <M 

*  ;  3 

Geometry:  y_  -  yV_./s  .  V_  -  200  m/s  (5) 

. . mm  m 

(the  hat  over  the  time  function  represents  its  Laplace  transform). 

aot2 

Input :  The  target  position  yTarg  -  yQ  ♦  vQt  ♦  — j—  (6) 

yn  €  I0.200]m  ,  vrt  €  (0,10]m/s 

°  2  °  (7) 

aQ  €  [0,2]m/s 

These  give  the  range  of  target  position  (initial),  velocity  and  acceleration 
values,  constants  in  a  specific  run. 


Feedback  equations: 

X  -  arc  tan  £*■  .  Ay  -  yTarg  -  ym  .  Ak  -  R  -  Vmt  cos  y  (8) 

R  (range)  6  [1500,3000]m  (9) 

N  (Fig.  2)  €  [2,M  (10) 

Initial  conditions: 

y(0)  -  y (0)  -  6(0)  -  0(0)  -  o(0)  -  6(0)  -  0  .  (11) 


Z 


j 


-7- 

3.  APPLICATION  Of  NONLINEAR  SYNTHESIS  TECHNIQUE 

3.1.  The  Set  of  Acceptable  Outputs  V  .  It  wes  emphasized  in  Section  I, 
that  It  is  imperative  that  a  set  of  acceptable  plant  outputs  V  ,  be  specified, 
In  order  to  obtain  the  tti  equivalent  set  P  .  The  latter  set  Is  guaranteed 
to  be  the  correct  equivalent  of  the  nonlinear  plant  set  (If  ,  only  for  the  set 
y  (although  it  may  and  usually  is  so,  for  a  larger  set  than  y  ).  Hence,  it 
is  essential  that  V  contain  the  true  typical  nonlinear  plant  outputs  which 
are  desired  under  actual  working  conditions,  including  the  extreme  cases. 

Here,  the  nonlinear  synthesis  technique  is  to  be  applied  to  the  portion 
labelled  A/P  (autopilot),  in  Figure  2  and  the  nonlinear  plant  (  w  of 
Section  l)  is  the  portion  from  6  to  y  .  So  the  set  V  refers  to  the  set 
of  desirable  y  that  exist  over  the  range  of  conditions  described  by 
Equations  (1-10),  with  particular  attention  paid  to  the  parameter  uncertainty 
spread  given  in  Equations  (A, 7, 9, 10). 

To  achieve  such  a  y  in  the  present  case,  the  -.utopi  lot  in  Figure  2  was 

*  * 

modelled  as  a  second  order  system  with  transfer  function  (for  y/A  ,  where 
the  hat  over  a  time  function  represents  its  Laplace  transform). 

:  Nui2 

f  *  -= - - - 2  '  N  of  E9-  HO)  (12) 

A  s  +  2^uj  s  +  a) 
n  n 

u>n  €  (1,4]  ;  t  =  .7  .  (13) 

This  second-order  simulation  of  the  autopilot  replaced  the  autopilot  in 
Figure  2.  Several  hundred  runs  were  made  of  this  simplified  Figure  2  system. 
The  fictitious  y(t)  values  so  obtained,  i.e.  the  output  of  this  second  order 
simulation,  provided  the  acceptable  set  V  .  Each  y  €  V  was  then  used  as 
a  known  signal  in  Equations  (la-c),  to  solve  for  0 ( t )  ,  <5(t)  ,  giving  an 
tti  equivalent  P(s)  . 

3-2.  The  LTI  Equivalent  Set 

Several  hundred  random  combinations  were  used  of  the  10  parameters  of 
Equations  (4,7,9,10,13)  in  the  above  manner.  For  each  combination,  the  time 
functions  6(t)  ,  0 ( t )  were  obtained  to  give  the  tti  plant  equivalent. 

Note  that  0 ( t )  rather  than  y(t)  is  used  at  this  stage  for  the  nonlinear 
plant  output,  because  the  feedback  is  taken  from  0  .  A  representative 
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nuefeer  of  these  time  function  pelrs  Is  shown  In  Figures  fce-q.  In  these 
end  later  figures: 


XN  -  N  (Eq.  10) 

9 

AT  -  aQ  (Eq.  7)  , 

2 

• 

< 

O 

m 

.0 

9 

XMO  -  yQ  (Eq.  7)  , 

• 

R  -  R  (Eq.  9) 

9 

WH  -  «n  (Eq.  13)  , 

0M 

CLA  -  Cu  (Eq.  3a) 

9 

CU>  “  CL6  (Eq*  3b) 

• 

CMA  -  C  (Eq.  3c) 

ma 

9 

CH0  “  Cm6  (Eq*  3d) 

• 

The  ttl  equivalent  P(jw) 

P(jw)  -  e(ju*)/6(jw) 

• 

05) 

3*3-  The  Frequency-Response 

Calculation  Problem 

At  first,  attempts  were  made 

to 

find  0(ju>)  ,  6(jw) 

separately  and 

perform  the  division  in  (15).  Adaptive  numerical  Integration,  of 

«• 

J  f(t)  dt  was  made,  in  the  sense  that  the  time  increment  for 

numerical  integration  was  shortened  as  u  was  increased.  For' many  cases, 
good  results  were  obtained.  For  many  others,  the  results  were  quite  poor  In 
the  high-frequency  range.  Other  computer  library  techniques, 
including  "Fast  Fourier  Transform",  were  tried  and  found  to  be  highly  in¬ 
adequate  In  the  high-frequency  range.  These  routines  are  adequate  for 
communication-oriented  problems,  where  good  accuracy  over  the  signal  band¬ 
width,  is  usually  sufficient.  However,  in  a  feedback-c.iented  problem  such 
as  the  present  one,  with  significant  parameter  uncertainty,  good  accuracy  Is 
needed  for  a  very  much  larger  frequency  range.  It  is  especially  Important 
to  have  clean  results  for  several  octaves  In  the  high-frequency  asymptotic 
region.  Thus,  given  a  signal  x(t)  and  Its  transform  X(s)  ,  as  s  ■*  •  , 
X(s)  ■*  k/s  p  ,  where  e^  Is  the  excess  of  X(s)  poles  over  zeros.  In  feed¬ 
back  problems  of  the  present  type.  It  is  essential  to  have  good  results  for 
several  octaves  in  this  range. 

In  our  problem,  It  turned  out  eventually  that  P(Jw)  ,  the  ratio  in 
05) ,  was  relatively  simple.  However,  each  of  6,6  was  very  complex,  with 
the  complexity  cancelling  out  In  the  division.  After  numerical  Integration 
was  found  Inadequate  at  high  frequencies,  It  was  thought  to  use  a  separate 
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representat ion  for  6(t)  ,  6(t)  ,  for  small  t  ,  since  the  latter  corresponds 
to  large  w  values  in  the  signal  transforms.  The  Taylor  series  expansion 
was  used.  A  computer  program  was  prepared  which  used  a  large  number  of  time 
samples  and  an  optimization  routine,  to  find  the  best  Taylor  series  fit  for 
some  range  t  €  [0,T/m]  ,  m  a  fraction  of  the  run  time.  Numerical  Integration 
was  also  used.  If  there  was  an  overlapping  frequency  range  where  the  two 
techniques  gave  very  similar  values,  then  the  results  were  considered  satis¬ 
factory.  The  numerical  integration  results  were  then  used  in  the  low  frequency 
range,  and  the  Taylor  series  results  in  the  high-frequency  range.  This 
approach  was  helpful  for  a  significant  number  of  runs,  but  there  still 
remained  a  substantial  number  of  cases  where  it  too  failed.  Figures  5a,b  are 
illustrative  of  the  | P ( j u>)  |  of  Equation  (15)  obtained  in  such  cases. 

Finally,  a  new  technique,  original  (to  the  best  of  our  knowledge),  was 
developed  in  which  P(jw)  was  found  directly,  without  first  finding  0(jw)  . 

6(j o>)  of  Equation  ( 1 5) .  In  this  way,  there  is  avoided  the  problem  of  finding 
the  very  complex  e(ju>)  ,  6(jw)  .  This  technique  is  described  in  Appendix  2, 
and  is  a  valuable  byproduct  of  the  research  effort.  Some  results  are  shown 
in  Figures  6a-f.  To  check  the  results,  numerous  runs  were  made  in  which  the 
<S  (t)  signal  was  the  input  to  Ui  network  with  the  transfer  function  P(jm) 
gotten  by  the  above.  Numerical  inverse  transform  integration  was  used  to  find 
the  resulting  output.  Representative  results  are  shown  in  Figures  7a-c.  In 
each  of  these  figures,  the  first  part  gives  the  8(t)  ,  6(t)  signals  obtained 
from  the  original  simulation  described  in  3-1*  The  second  part  shows  |P(jw)| 
obtained  by  the  new  method.  Arg  P(jw)  was  obtained  but  is  not  shown.  The 
third  part  gives  L  P(jw)6(jw)  *  6*(t)  ,  which  should  be  checked  against 
0(t)  of  the  first  figure.  It  is  seen  that  excellent  agreement  was  obtained. 

3-1*.  Review  of  LTI  Synthesis  Technique 

Once  C  ,  the  fcti  plant  equivalent  of  the  nonlinear  plant  set  W  has 
been  obtained,  the  problem  becomes  one  of  designing  a  tti  feedback  system 
to  cope  with  uncertainty  of  a  tti  plant.  A  detailed  synthesis  technique 
for  this  purpose  is  described  in  Reference  3,  which  is  briefly  reviewed  here. 
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Figure  8  is  used,  in  which  the  closed-loop  transfer  function  is 


T(s)  -  F(s) 


It  is  assumed  that  the  compensation  networks  F  ,  G  , 


' ’  '  1+GP(s)  •  '  - - -  -  - - -  ’  ’ 

whose  power  level  can  be  very  low  (as  the  plant  contains  the  power  elements), 

can  be  constructed  with  negligible  uncertainty  in  their  transfer  functions. 

Hence,  due  to  the  uncertairity  in  P  , 


AAnT  -  A  In  -  Atn  , 

Atn|T(Jw))  -  Atnl-j-^J^I  . 


l-GP 


There  are  given  u-domain  specifications  on  |T(ju) |  ,  as  in  Figure  9.  Given 
that  the  maximum  allowed  change  A£n|T(j<<>)|  (  6j  db  for  example  at  in 

Figure  9  ,  what  are  the  resulting  constraints  on  L(jw)  7  It  is  convenient 
to  pick  a  "nominal"  plant  Pq(s)  ,  and  derive  the  bounds  on  the  resulting 
"nominal"  loop  function  Lq  -  PgG  .  These  bounds  can  be  found  by  means  of  a 
digital  computer,  but  it  is  very  useful  for  insight  to  see  it  done  on  the 
Nichols  chart  (logarithmic  complex  plane  with  abcissa  in  degrees,  ordinate 
in  decibels  *  20  log^) .  The  procedure  is  illustrated  for  the  case 


P(s)  *  TTifsy  ;  1  1«a<10  . 


At  u»2  rps  ,  P(j2)  lies  within  the  boundaries  given  by  ABCD  in  Figure  10. 
Since  AnL  ■  AnG  +  AnP  ,  the  pattern  outlined  by  ABCD  may  be  translated,  but 
not  rotated,  on  the  Nichols  chart,  the  amount  of  translation  being  given  by 
the  value  of  AnG(j2)  .  For  example,  if  a  trial  design  of  L ( j 2 )  corresponds 
to  the  template  P(j2)  at  A'B'CD'  in  Figure  10,  then 


I«(i2)ldb  -  IUJ2)|db  -  |PtJ2)|d(> 

-  (-2.0)  -  (-13.0)  »  11.0  db 
Arg  G(j2)  -  Arg  L(j2)  -  Arg  P(J2) 

-  (-60*)  -  (-153. *’)  -  93- k#  • 


(19a) 


(19b) 


3.k.2.  Bounds  on  L(lu)  in  the  Nichols  chart 

The  templates  of  P(jw)  are  manipulated  to  find  the  position  of  L(jw) 
which  results  In  the  specifications  of  Figure  9  on  £n|T(ju)|  being  satisfied. 
Taking  the  u*2  template,  one  tries,  for  example,  ;-ositioning  it,  as  shown 


R(s) 

r(tT 


FIGURE  8 
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d(t) 


LTI  two-degree-of “freedom  structure. 


FIGURE  9  Bounds  on  |t<1w)!» 


in  Figure  10,  at  A'B'C'O'.  Contours  of  constant  £n|L/(l+L)|  are  available 
on  the  Nichols  chart.  Using  these  contours,  it  is  seen  that  the  maximum 
change  in  £n|L/(1+L)|  ,  which  from  (17)  is  the  maximum  change  in  fcn|T|  is, 
in  this  case,  very  closely  (-0.49)  -  (~5.7)  ■  5.2  db  ,  the  maximum  being  at 
point  C',  the  minimum  at  point  A*.  Suppose  that  the  specifications  tolerate 
a  change  of  6.5  db  at  u>  -  2  ,  so  the  above  trial  position  of  |  L  ( j  2 ) j  is  in 
this  case  more  than  satisfactory.  The  template  is  towered  on  the  Nichols 
chart  to  A"B"C"0",  where  the  extreme  value  of  £n|L/(l+L)j  are  at  C" 

(-0.7  db).  A"  (-7-2  db).  Thus  if  Arg  Lft( j2)  -  -60°  ,  then  -4.2  db  is  the 
smallest  magnitude  of  ( j 2 )  which  satisfies  the  6.5  db  specification  for 
Aln|T|  .  Any  larger  magnitude  is  satisfactory  but  represents  over-design  at 
that  frequency.  The  manipulation  of  the  w *  2  template  is  repeated  along  a 
new  vertical  line,  and  a  corresponding  new  minimum  of  |L^(j2)|  found. 
Sufficient  points  are  obtained  in  this  manner  to  permit  drawing  a  continuous 
curve  of  the  bound  on  ( J 2 )  ,  as  shown  in  Figure  10.  The  above  is  repeated 
at  other  frequencies,  resulting  in  a  family  of  boundaries  on  ( j 2 )  . 

3.4.3.  Nature  of  the  bounds  on  L ( j u>) 

A  typical  set  of  bounds  is  shown  in  Figure  11.  The  bounds  tend  to  move 
down  in  the  Nichols  chart  (become  less  onerous),  obviously  because  as  w 
increases,  greater  change  in  ! T ( j u>)  j  is  permitted,  as  in  Figure  9-  It  is 
in  fact  essential  that  at  large  enough  u>  ,  the  uncertainty  in  | T ( j cu)  | 

(i.e.,  the  bounds  on  jT(jui)  |  )  be  greater  than  the  uncertainty  in  P(ju>)  , 
because  the  net  sensitivity  reduction  is  always  zero  in  any  practical  system 
as  was  long  ago  (4)  shown  by  Bode, 

00  ao 

J  »n jsl(ju>)  |dui  *  -  J  £n|  1  +  L(ju>)  |du>  =  0  (20) 

0  0 

where  sj  »  sensitivity  function. 

in  the  above  example,  as  w  -*■  «  ,  P  -*■  ka/s*  ,  so  A£n|P|  -*  AJtn(ka)  - 

40  db  .  Note  in  Fiqure  9  that  the  permitted  A£n|T(jw)|  >>  40  db  for 

u>  >  50  .  Such  large  tolerances  on  J T ( J u>)  |  at  large  u>  are  tolerable  because 

|  T  ( j  t*>)  |  is  negligible  at  large  w  ,  e.g.,  if  |  P  ( j  u>)  |  can  change  at  most  by 

40  db  at  large  w  but  |T(jw)|  changes  by  52  db,  who  cares  if  this  52  db 

change  is  from  I T I  .  *  10  ^  to  I T I  to  400  *  10  In  return,  one  can 

*  1  'mm  '  'max 


DEGREES 


Typical  bounds  on  Lo(ju>)  at  low  and  medium  frequencies 


concentrate  the  sensitivity  reduction  over  the  bandwidth  of  T(jui)  .  Thus, 
although  |P(iw)|  in  this  region  vales  by  say  40  db,  jT(jw) |  may  be 
controlled  to  vary  by  only  4  db,  or  0.04  db  if  desired. 

3.4.4.  Universal  high-frequency  boundary 

As  noted,  in  the  high-frequency  range  Atn|T(jw)|  must  realistically 

be  allowed  to  be  >>  Ain|P(ju))|  ,  and  this  is  reflected  In  the  bounds  on 

Lg(jw)  tending  to  a  very  narrow  pencil.  In  Figure  12,  Is  drawn  for  the 

case  At.nl  ■  Atnk  -  20  db  ,  AtnjT{ju>)|  -  Atn  ( L/  ( 1+1)  )  *  23  db  at  u>v  . 

However,  the  resulting  peak  value  of  |l/(l+L)|  is  23  db  ■  14.1  arithmetic 

at  4  "  ’  Indicating  a  highly  under 'damped  pole  pair  at  the  corresponding 

frequency  with  damping  ratio  f.  "  0.034  ,  when  4  *  k  ^  .  This  tremendous 

peaking  does  not  appear  In  the  system  response  to  the  command  inputs  R  , 

because  it  is  filtered  out  by  the  pre-filter  F  in  Figure  8.  But  the  system 

response  to  a  disturbance  0  is  given  by  »  C/0  ■  (l+l)  *  .  Disturbance 

attenuation  generates  its  own  requirements  on  l  ,  which  may  lead  to  more 

stringent  bounds  on  l  than  those  due  to  T  (  J  <»\)  .  The  final  contours  used 

in  the  design  1 3 1  must  be  the  most  stringent  composite  of  the  two.  However, 

even  if  D  is  very  small,  It  is  usually  certain  that  a  peak  |l^|  of  l4.1 

Is  intolerable.  It  is  reasonable  to  add  a  requ i i oment  | |  (y  some 

constant,  for  all  u>  and  over  the  whole  range  of  P  parameter  values. 

The  resulting  constraining  contours  denoted  by  are  shown  in  Figure  12 

for  the  case  Atnk  -  20  db  ,  and  for  y  •  2.3,  3.5.  b  db  (all  these  contours 

are  symmetrical  with  respect  to  the  verical  line  Arg  l  »  -l80°  on  the 

Nichols'  chert).  If  y-5  db  is  used,  then  B(w  )  indicates  the  composite 

contour  shown  in  Figure  12.  For  ui  '  uiy  ,  |  A  T  ( |  u> )  |  increase  while  y 

remains  the  same,  so  that  sooner  or  later  there  Is  i cached  a  frequency 

w  3  B(w)  is  a  fixed  boundary  B^  ,  effective  V  m  -  m 

y  h  y 

3.4.5.  The  optimum  l(jdi) 

It  has  been  shown  |5|  a  realistic  definition  of  optimum  in  the  t.ti 
system  is  the  minimization  of  k  ,  defined  by  1 1  in  l(s)  -  ks  e  ,  where  e 

is  the  excess  of  poles  over  zeros  assigned  to  L(s)  . 

It  has  been  proven  (51  that  the  optimum  L  lies  on  its  boundary  Bj 
at  each  «j  and  that  such  an  optimum  exists  and  is  unique.  Most  Important 
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for  the  present  purpose,  is  that  in  significant  plant  Ignorance  problems 
the  ideal  optimal  L  has  the  properties  shown  in  Figure  13,  i.e.,  over  a 
significant  range  it  follows  along  UV  to  the  point  J  at  which  it 

abruptly  jumps  to  infinity  along  WW'W"  and  returns  on  the  vertical  line  Y Z, 
whose  phase  is  (-90°)*e  .  Such  an  ideal  L(jw)  is,  of  course,  impractical. 
A  practical  suboptimum  L  is  shown  in  Figure  1 3 • 

Some  results  of  a  numerical  design  example  are  shown  in  Figure  14. 

They  were  derived  for  the  following  problem. 


Numerical  example. 
Plant : 

Plant  ignorance: 

Performance  Specificati 

Disturbance  response: 


(Figure  14) 

P  -  PlP2 

P,  =  kj/s  ,  /2  s  *  10/2 

P2  =  k2/s  ,  /2  f  k2  s  \0/2 

on:  Shown  in  Figure  9  were  originally  derived 
from  time  domain  bounds  15). 

y  s  2.0  db  . 


The  derivation  of  a  rational  Lq(s)  which  satisfies  the  bounds  B(u) 
on  Lq ( j  to)  ,  is  somewhat  of  an  art.  For  a  given  skill  in  the  art,  there  is 
trade-off  between  design  complexity  and  Lq ( j u>)  bandwidth  reduction.  Some 
computer  proqrams  have  been  prepared  for  automatic  derivation  of  a  pinctical 
optimum  L^tim)  ,  but  complicated  B(w)  patterns  require  human  interaction. 

3.4.6.  Des i gn  of  pref i I  ter  F (s) 

The  above  only  guarantees  that  A  log  |T(ju>)  |  s-  the  relative  change 
allowed  by  the  tolerances  in  Figure  9.  For  example  at  Wj  ,  say  the  allowed 
^max  *  db  ,  and  allowed  l^lmjn  *  "22  db  ,  i.e.  a  maximum  change  of 

14.4  db.  A  proper  design  of  Lq(s)  in  the  manner  described  above  only 
guarantees  that  the  change  of  j L/ l+L |  *  14 . 4  db  .  But  it  is  possible  for 
|L(jhi|)/l  ♦  L  ( )  i*>  | )  |  to  actually  be,  for  example,  in  the  range  1.2  db  to 
-8  db.  The  function  of  the  prefilter  F(jw)  in  Figure  8,  is  to  achieve  the 
shift  needed.  In  this  example  any  |F(ju>j)|  €  [-8.8,-l4)  db  is  satisfactory, 
an  allowance  of  5.2  db  because  LqQwj)  was  overdesigned  by  5*2  db. 


L  L  (degree) 


Practice!  euboptimum  l 


FIGURE  1J  Optimum  and  "praotiraj  opt  imum"  (,. 


L L»  (DEGREE) 


FIGURE  14  A  dos  iqn  oxamplo. 


I  L(  (db) 
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We  now  turn  to  the  application  of  the  methods  of  Section  3-^  to  our 
present  problem. 

3.5*  Design  of  Autopilot 

Section  3-3  described  how  the  frequency  response  P ( j u>)  of  the  Iti 
equivalent  plant  was  obtained.  Following  3.^.1 .  the  next  step  is  to  find 
Tp(u>)  templates  of  P(jui)  ,  for  a  reasonable  number  of  discrete  w  values. 
Tp(b>|)  is  simply  the  set  I  P(  ju^ )  1  .  A  number  of  these  are  shown  In 
Figure  15*  in  which  the  case  marked  by  the  large  star  is  the  chosen  nominal 
plant  function  Pg(ja>)  .  At  very  small  m  ,  Tp(w)  is  a  vertical  line  because 
then  P(jw)  ►  a  constant  of  /eio  phase  but  a  different  constant  for  different 
parameter  combinations.  At  large  u  ,  each  P(s)  >  kp/s  in  the  case,  with 
kp  a  function  of  the  parameter.  Hence  Tp(m)  is  again  a  vertical  line. 

Note  the  relatively  small  area  of  Tp(u>)  ,  indicating  rather  small  un¬ 
certainty.  Compare  these  Tp(to)  with  those  in  the  aircraft  flight  control 
problem  of  Appendix  1  (see  Figures  6a-c,  Appendix  1).  In  the  latter  some 
P  €  P  are  unstable  (e.g.  cases  e,f  in  Figure  5  there).  This  is  also  seen 
from  Appendix  1,  Figure  6,  where  there  are  two  verical  line  sub-templates 
360°  apart,  due  to  t£e  two  groups  of  Iti  P(jw)  ,  one  stable  and  the  other 
with  two  right  half-plane  poles.  The  uncertainty  and  nonlinearity  in  this 
missile  problem  are  considerably  less  than  in  the  aircraft  problem  of 
Append i x  1 . 

The  bounds  on  |  T  ( j  *■*)  |  (recall  3.*t,  Figure  9)  are  needed.  These  are 
determined  by  the  model  T ( s )  in  Equations  (12,13),  used  for  the  autopilot. 
The  resulting  permissible  spread  in  |T(ju>)|  ,  due  to  allowed  to  range 

€  [1*M  ,  is  shown  in  Figure  16.  The  procedure  described  in  3.^.1  (Figure  10, 
etc.)  is  next  followed  to  find  the  bounds  on  LqOm)  .  These  are  stiown  in 
Figure  17*  The  bounds  are  very  moderate  ones,  due  to  the  relatively  small 
Tp(u>)  sizes  (small  effective  uncertainty  in  P(jw)  ),  and  generous  | T ( j u>)  | 
to  I  erances.  In  view  of  these  small  demands,  a  very  generous  Lq ( j <o)  was 
chosen  which  satisfied  the  bounds  B(w)  ,  with  much  to  spare  —  shown  in 
Figures  17,18. 
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Oesiflo  of  Prefilter.  The  design  of  F(s)  of  Figure  8,  to  satisfy  the 
bounds  on  J T < ju>)  J  of  Figure  9,  has  been  described  in  3-4.6.  Our  problem 
is  more  complicated  because  the  feedback  loop  is  around  e  ,  while  the 
tolerances  are  on  the  y  response  in  Figure  2.  It  is  therefore  necessary 
to  derive  0/y  which  multiplies  the  bounds  on  y  ,  to  give  the  resulting 
bounds  on  8  .  A  serious  difficulty  would  arise  if  the  uncertainty  (spread) 

*  '  It  A 

in  0/y  was  of  the  same  order  of  magnitude  as  that  of  P(jw)  »  8/6  .  It  is 
essential  that  the  former  be  somewhat  less  than  the  latter.  Fortunately,  in 
this  example  there  is  little  uncertainty  in  0/y  .  Numerous  runs  (as  des¬ 
cribed  in  3*2),  were  made,  giving  0(t)  ,  y(t)  .  The  new  technique  of 

A  A 

Appendix  2  was  used  to  find  0/y  .  Representative  samples  are  shown  in 

Figures  I9a-f,  and  it  is  seen  that  there  is  little  uncertainty  in  0/y  .  In 

each  of  Figures  19a-f,  the  first  part  gives  the  time  functions,  and  the 
second  the  frequency  response  1 8( jw) /y (jw)  |  .  The  new  bounds  on  |T(jw) |  , 

due  to  0  being  considered  as  the  autopilot  output,  are  given  in  Table  2. 

The  old  bounds,  due  to  y  as  output,  are  given  in  Table  1.  The  bounds  on 
| F ( j u»)  |  were  then  obtained,  as  explained  in  3-4.6,  resulting  in 

F(s) "  (s-+T:5y(s'+M(s  +  8) 


3.6.  Design  Simulation  —  Local 

It  is  important  to  recall  that  the  nonlinear  design  technique  has  been 
applied  only  to  the  autopilot  in  Figure  2  —  not  to  the  entire  closed- loop 
system  of  Figure  2.  The  proper  verification  of  the  design  is  then  only  for 
this  autopilot  portion.  To  effect  this  verification,  the  simulation  of  3.2 
was  used  to  derive  X  which  became  the  input  driving  our  above-designed 
autopilot.  In  the  latter,  there  were  used  aerodynamic  equations  (1-3)  and 
the  G(s)  ,  F(s)  of  3-5.  Thus,  only  the  autopilot  as  a  self-contained 
closed-loop  system  was  checked  at  this  point  —  not  as  a  part  of  the  larger 
closed-loop  structure  of  Figure  12.  The  results  of  the  simulation  for  a 
representative  number  of  runs,  including  the  extreme  (largest  0  values), 
are  shown  in  Figures  20a-e,  for  Q ( t )  and  6(t)  .  Each  figure  has  three 
parts.  The  first  repeats  the  result  obtained  in  the  simulation  of  3*2  with 
the  second-order  T(s)  model  of  Equations  (12,13)  (i.e.  same  as  in 
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Figures  (Aa-q)).  The  second  part  uses  the  same  input  A  (to  the  autopilot) 
as  in  the  first  part  and  the  outer  loop  is  not  closed.  It  is  only  a  test  of 
the  autopilot  and  Is  therefore  denoted  as  "local  feedback  simulation".  In 
the  third  part  the  outer  loop  is  closed,  so  the  simulation  is  that  of  a 
"practical"  system. 

It  is  seen  that  in  general  the  6(t)  ,  6 ( t)  values  in  the  last  two 
parts  are  larger  in  magnitude,  but  considerably  smoother  than  in  the  first. 

But  in  considering  the  actual  tracking,  the  outer  loop  was  not  found 
satisfactory.  This  is  next  treated. 

3.7*  Outer  Loop  Modification  and  Global  Simulation 

The  above  autopilot  design  was  next  Imbedded  Into  the  over-all  structure 
of  Figure  2,  and  the  tracking  examined  (  y^  and  Ym  ) •  In  a  sense  it  is  an 
unfair  comparison  of  the  nonlinear  design  technique,  because  the  latter  was 
executed  only  for  the  autopilot  part.  In  any  case,  when  effected,  the  poorest 
results  were  the  three  shown  in  Figure  21,  in  which  yT(t)  represents  the 
target  trajectory  and  ym(t)  .  the  missile  trajectory.  The  run  is  stopped 
precisely  when  R  ■  Vmt  cos  y  ,  i.e.  when  Ax  ■  0  .  A  perfect  hit  requires 
Yj  “  ym  precisely  at  this  end-point.  In  Figure  21a,  there  is  an  error  of 
6  meters,  and  much  smaller  errors  in  21b, c.  Also,  the  ym(t)  trajectories 
In  Figures  21a, b  are  much  too  oscillatory. 

The  reson  for  the  above  is  due  to  the  fact  that  the  design  was,  up  to 
this  point,  completely  a  local  one  —  of  the  autopilot  alone,  to  achieve  a 
specified  response  range  for  the  autopilot.  No  overall  global  closed- loop 
response  was  considered  at  all.  The  logical  procedure  at  this  point  would 
be  to  find  the  tti  equivalent  of  the  time-varying  portion  of  the  structure 
(i.e.  of  the  A  •  tan*1  Ay/Ax  portion)  and  apply  the  design  technique  to 
the  outer  loop.  The  nonlinear  design  technique  applies  also  to  linear  and 
nonlinear  time-varying  uncertain  plants  —  see  Reference  1.  However,  there 
was  insufficient  time  left  for  effecting  such  a  design.  Therefore  an  ad-hoc 
modification  of  the  prefilter  F(s)  was  made,  with  the  new 
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This  gave  much  better  results,  shown  in  Figures  22a-t,  whose  inter¬ 
pretation  is  the  same  as  of  Figure  21.  The  star  (explosion)  indicates  an 
"error"  <  2  meters,  Interpreted  as  a  "hit".  Of  ail  the  runs  made,  there 
was  only  one  case  where  the  error  ft*  3  meters,  in  Figure  22s. 

To  complete  the  picture,  a  large  number  of  runs  of  the  variables 
6,o,  accel  ■  Vy/g  are  shown  in  Figures  23a-p. 
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AI'I'FNHIX  1 

FLIGHT  CONTROL  OKS  ION  RASFO  ON  NONLINI.AK 
MODF.I.  WITH  HNCKRTAIN  PARAMKTF.RS 

I  *  * 

Isaac  Horowitz  ,  Ho*  is  Golubev  ,  Tzvi  Kopclman 

Hop,**  t  men  I  of  Applied  Mat  l»  om.it  i  os 
The  Wotzm.inn  Institute  of  Selenee 
Rehovot ,  Israel 

Aba  t tar t 

A  simplified  but  s  i  gn  i  f  i  e.»nt  l  v  noitiine.i*  model  1  •.  n:ed  to*  the  she*  t  petiod 
longitudinal  flight  control  problem,  with  cA  taken  as  the  output.  Uneeitaintv 
la  included  in  the  model  by  allowing  for  a  large  range  in  velocity  and  air-dcusitv, 
without  any  provision  lor  their  measurement.  The  output  o*(tl  is  to  lie  within 
specified  bounds  in  response  to  a  range  of  step  commands.  A  recent  exact  design 
technique  for  uncertain  nonlinear  systems  is  used.  In  this  technique,  the  eon 
linear  plant  set  Is  replaced  bv  a  linear  time  Invaiiant  plant  set,  which  is 
precisely  equivalent  to  the  nonlinear  set  with  respect  to  the  problem  toleiances. 
The  design  execution  then  involves  frequency  response  concepts  and  techniques. 

For  a  large  problem  class,  the  resulting  design  is  guaranteed  to  solve  the 
problem,  for  the  specific  class  of  inputs  posed.  The  design  was  simulated  with 
excellent  results. 

Index  Categories:  Guidance  and  Control;  Handling  Qualities,  Stability  and 
Cent  rol . 

t  Cohen  Professor  of  Applied  Mathematics 

*  Studenta 
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Nomenclature 


a  -  parameter  in  model  of  T(s) 

B(u>)  -  bounds  on  G(jw) 

C  —  mean  aerodynamic  chord,  m  . 

C^j(a)  various  aerodynamic  coefficients 

c*  —  output  variable  =  a  handling  qualities  criterion 


C*(s) 


c*(t) 


db  —  decibels  (20  log^) 
e  -  excess  of  poles  over  zeros  of  G(s) 

c 

F(s)-  prefilter  transfer  function 

2 

g  —  gravity  acceleration  9.8  m/sec 
G(s)  -  loop  compensation  function 

2 

Xy  —  moment  of  interia  in  pitch,  kgm-m 

j  -  /=T 

k  —  see  R 

lti  -  linear  time-invariant 
L(a)  »  G(s)P(s)  ,  loop  transfer  function 
-  Laplace  transform 
m  —  mass  of  vehicle,  kgm 
mp  —  minimum-phase 
nmp  -  non  minimum-phase 
p(t),  Pj  equivalent  lti  plant  function 

PC)  j£p(t> 

jP  set  { p>  or  (P(s)} 

q  ■  8  -  pitch  angular  velocity 
r(t)  system  command  input 
R(s)  j£r(t)  ”  9 


set  (r(t))  or  (R(s)} 
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comp  lex  variable 

2 

wing  surface  area,  m 

closed  loop  transfer  function 

set  {T(s) } 

template  of  P  -  P(jw) 

vertical  velocity,  m/sec 
horizontal  velocity  -  VQ  constant 
w  —  nonlinear  plant  function,  y  -  w(x) 

w  1  -  inverse  of  w  ,  x  -  w_1(y) 

set  (w) 

x,  Xj  input  to  nonlinear  plant 
X(s)  -  <£x(t) 

y. V^  acceptable  output  of  plant 
*<s)  J[y(t) 

^  set  fy} 
o  angle  of  attack 

5  elevator  deflection 


C  parameter  in  model  of  T(s) 

p  air  density,  kgm/m3 

u>  frequency,  rps 


S  — 
T(s) 

3- 

3- 

u  — 

V  — 


FLIGHT  CONTROL  DESIGN  BASED  ON  NONLINEAR 


MODEL  WITH  UNCERTAIN  PARAMETERS 


INTRODUCTION 

In  most  flight  control  design  techniques,  the  nonlinear  differential 
equations  are  linearized  about  a  trim  condition.  The  resulting  Incremental 
linear  time  Invariant  (lti)  model  with  its  fixed  aerodynamic  coefficients. 

Is  reasonably  valid  for  small  excursions  from  trim.  By  using  different  trim 
points,  there  are  obtained  different  sets  of  coefficient  values.  The  model 
is  then  taken  as  lti  ,  in  which  these  sets  are  treated  either 

(1)  as  uncertain  parameters  in  a  lti  design  technique  for  securing  specified 
performance  over  the  uncertain  parameter  set,  or 

(2)  as  known  functions  of  Mach  and  dynamic  pressure.  The  latter  are 
monitored,  leading  to  a  "Scheduling"  design,  or 

(3)  as  uncertain  parameters  to  he  identified  and  compensated  for  In  a  so-called 
"Adaptive"  design. 

A  combination  of  all  three  may  also  be  used.  However,  in  all  of  these 
approaches  the  model  used  is  It i  . 

There  have  recently  been  attempts  at  incorporating  the  non l inear i t tes  to 

acme  extent  into  the  system  model,  the  vehicle  for  doing  so  being  optimal  control 

theory.  However,  the  calculations  have  been  very  laborious,  requiring  oonslderabl 

approximations  (see  Ref.  I  and  its  references  for  discussion).  This  paper  takes 

2 

a  far  different  approach  and  is  based  on  a  recent  synthesis  technique  for 
feedback  aronnd  a  nonlinear  uncertain  plant,  the  latter  denoting  the  constrained 
part  of  the  system.  This  synthesis  technique  has  the  following  important 
properties. 


1.  For  a  large  class  of  practical  nonlinear  plants,  the  design  is  precise 
and  direct  with  no  approximations  even  for  highly  nonlinear  (even  nonlinear 
time-varying)  models  with  large  uncertainties  in  the  plant  parameters. 

2.  Design  execution  is  in  the  frequency  domain. 


REVIEW  OF  NONLINEAR  DESIGN  TECHNIQUE 

The  key  ingredient  is  the  extension  to  sets  of  an  idea  often  used  by 
engineers:  Given  a  specific  nonlinear  element  w  ,  it  is  usually  possible  to 
find  an  Iti  element  P  which  is  equivalent,  in  a  certain  sense,  to  w  for 
a  single  specific  output  y^(t)  (or  Input  x^(t)):  Let.  y^(t)  »  w(x^(t)) 
be  the  (assumed)  unique  output  of  w  due  to  input  x^(t)  .  Choose  p  so 
that  its  output  is  also  y^  when  its  input  is  x^(t)  .  A  simple  way  to  do 

(t)  -  Y1(s)  and  x  j  ( t )  »  X^s)  . 

Then  let  P(s)  *  Y^(s)/X^(s)  ,  giving  p(t)  '  ^  P(s)  the  equivalent 

of  w  t  but  only  for  the  special  case  when  the  input  is  x^  .  Hence,  the  two 
essentia]  conditions  are:  (1)  w  has  a  unique  inverse,  thus  excluding  hard 
saturation  etc.  which,  however,  could  be  replaced  by  very  low  gain  over  the 
appropriate  interval.  (2)  the  desired  plant  output  y(t)  and  the  resulting 
x  •  w  *(y)  ,  are  Laplace  transformable,  a  condition  difficult  to  violate. 


this  is  to  find  (Laplace  transform) 


l 


This  eqi valence  idea  can  be  extended  to  a  set  of  nonlinear  plants 
(w)  for  which  a  Iti  set  denoted  by  can  be  found  which  in  equivalent 

to  »  with  respect  to  y^(t)  .  Simply  find  p^  the  y^-equivalcnt  of  w^  , 
for  each  w^  €  and  let  ft-  {p^}  .  A  further  extension  is  to  find  a  Iti 
•«t  ^ which  is  eqivalent  to  JY with  respect  to  a  set  of  outputs  { Yj  > . 


f. 


repeat  the  previous  for  each  y.  €  Y  giving 


and  then 


L 


To  find 


e-  </?  >.  -  ^  has  10  elements  and  ^  has  20  ,  has  200  members. 

But  in  general,  both  X-  and  are  uncountable  and  so  is  . 

In  Fig.  1,  the  closed-loop  system  Is  to  behave  like  a  lti  one,  with 
transfer  function  T(s)  ,  in  response  to  any  command  input  in  a  set  =  {r) 
and  for  any  w  -  i.e.  y(t)  ^  Y(s)  =  T(s)R(s)  ,  where  R(s)  ■=  l  r(t)  . 

However,  due  to  the  uncertainty,  a  sot  of  acceptable  {T(s)}  -  ^  must  be 
specified,  insomuch  as  dynamic  invariant  response  for  all  w  in  >  is 
impossible.  The  sets  determine  the  set  of  desired  acceptable  system 

outputs  «  { y ( t)  }  ,  via  /,<  t)  »  Y(s)  =  T(s)R(s)  ,  or  vice  versa.  The 
nonlinear  set  is  then  replaced  by  the  lti  set  ^  which  is  equivalent 

to  %  with  respect  to  ^  .  We  say  that  iP  is  mp  (minimum-phase)  with 

respect  to  ^  ,  if  all  P(s)  €  ^  have  no  zeros  in  the  right  half-plane, and 
stable  if  they  all  have  no  such  poles.  In  our  specific  problem  (p  is  mp  , 


but  not  stable. 


We  now  have  a  pure  lti  feedback  problem:  Find  the  lti  compensation 
F(s)  ,  G(s)  needed  in  Fig.  1,  such  that  the  system  transfer  function 
T  •  Fl./O+L)  ,  l.  *  GP  ,  is  a  member  of  the  set  ,  no  matter  which  P  €  (P 

is  used.  A  frequency-response  design  technique  for  a  very  general  problem 


class  ,  is  available  for  minimum-phase 


e.  e 


may  have  unstable  elements. 


If  this  lti  problem  is  solvable,  then  under  quite  general  conditions  the 

solution  (i.e.,  the  F(s)  ,  G(s)  compensation  pair  used)  is  also  precisely 

2 

valid  for  the  original  nonlinear  problem.  The  details  and  proof  involve 
functional  analysis  techniques.  However,  design  execution  involves  frequency- 


response  concepts,  as  next  seen. 
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PROBLEM  STATEMENT  AND  EXECUTION 

Nonlinear  Plant.  The  nonlinear  plant  is  described  by  Kqs.  1-3  (see 
Nomenclature  Section). 

pVos 

u  »  qVQ  +  g  cos  0  -  -gj-  (CNa(a)  +  CN6(a)6]  (1) 

pV  SC 

0  -  4  »  “21^-  f Cm6(c)«  +  (^a(a)  +  2 ^  Cmq(a)q)  .  a  =  tan"  ^  .  (2,3) 

c*  “  12.40  +  “  (u  -  vJ  +  60)  (4) 


9.8 


is  the  output  to  be  controlled,  i.e.  to  be  i 


n 


The  numbers  used  are 


5.6 


2  2 
1^  *■  207,000  kg  m  ,  m  =  17,600  kg  ,  C  =  4.89  meters  ,  S  -  49.2  m  .  The 

Cjj(a)  are  nonlinear  functions  of  a  ,  see  Fig.  2«  Since  a  (Fig.  9  ) 

ranges  in  10,35°],  there  is  strong  nonlinear  operation.  The  horizontal 

velocity  v  was  taken  as  Vq  fixed,  which  is  incorrect  for  some  low-velocity 

cases,  but  the  objective  here  is  to  demonstrate  the  validity  of  the  design 

technique  in  a  strongly  nonlinear  situation,  which  is  achieved  sufficiently 

by  means  of  the  nonlinear  C^(a)  . 

The  bounds*1  on  the  acceptable  c*(t)  in  response  to  a  unit  step  command 
are  included  in  Fig.  3,  which  also  includes  design  simulation  results.  The 
set  of  command  inputs  consists  of  steps  1  to  5  in  magnitude.  Parameter 

uncertainty  is  due  to  0  ranging  in  [.3,1.22]  and  Vq  in  [75,206]  .  Initial 
conditions  are  w(0)  ■  q(0)  *  q(0)  »  0  ,  a(0)  >*  0(0)  ,  giving  initial  values 
for  6  (as  well  as  u  ,  0  ,  a)  which  is  subs trac ted  out  so  that  the  change 
in  6  is  used  to  find  the  lti  equivalent  set  .  The  detailed  steps  in  the 

design  implementation  are  next  presented,  with  comments  postponed  to  the  end. 


I 


k  *  The  lti  set 


y  .  Let  ^c*(t)  -  C  *  ( s  )  -  T(s)R(s)  , 


R(s>  -  k/s  with  k  €  [1,5]  and  T(s)  €  ^  ,  derived  from  the  bounds  in 
Fig.  3.  A  simple  means  for  generating  ,  taken  from  Ref.  5,  Is  to  let 
T(s)  ■  a^(s  +  2.9)/2.9(s^  4  2Cas  4-  a4-)  ,  C  ranging  in  [3.7,  1.5]  ,  a  in 

(3.14,  7.6]  ,  giving  the  bounds  on  |  T  (  (»i)  |  in  Fig.  4.  Such  hounds  suffice' 

for  mp  (and  obviously  stable)  T(s)  .  Any  T(s)  ,  R(s)  pair  thus  generates 

* 

an  acceptable  c  (t)  .  A  computer  program  solved  Fqs .  (1-4)  backwards  for 
6(t)  and  then  chocked  the  result  hv  solving  F.<|s .  (1,4)  forwards  for  c  (t) 

from  A(t)  .  The  program  was  considered  adequate  only  when  excellent  agreement 

* 

was  obtained  over  the  entire  range  considered.  C  (s)  was  apriori  available  and 

A(l.ii)  ”  was  obtained  by  numerical  i nteurat ion.  As  c*  is  a  short -per iod 

3 u' 

criterion  and  for  all  the  acceptable  cases  has  definitely  reached  steady-state 
in  4  seconds  (see  Fig.  1),  Fqs.  (1,4)  were  solved  only  for  t  £  [0,4  seconds], 
and  the  constant  5(4)  was  used  for  t  '  4  . 

Loci  of  six  P ( J to)  are  shown  in  Fig.  5,  two  of  them  (e,f)  unstable 
with  a  pair  of  right  half-plane  poles,  which  are  zeros  of  A(s) .  The  set 
includes  a  large  number  of  such  unstable  lti  1’ts)  ,  which  the  lti  design 
technique''  can  easily  handle. 


Plant  templates.  At  any  cu  say  io  -  w ^  ,  the  set  { P ( J»o> )  ,  P  £  jg 
consists  of  a  region  in  the  logarithmic  complex  plant  (Nichols  chart)  denoted 
as  the  u<j-plant  template  ^Jp(tOj)  •  A  number  of  ^(u>)  are  shown  in 
Figs.  6a-e.  At  very  small  w  there  are  two  almost  constant  angle  sub- templates 
360°  apart.  This  Is  due  the  presence  of  both  stable  and  unstable  P  £  (P  . 


and  the  fact  that 


Arg  P  near  u>  -  0  is  either  *  0  or  nn/2  for 


some  Integer  n  .  As  u>  increases,  the  two  groups  merge  together  and 
approach  a  single  vertical  line  at  large  w  well  beyond  the  plant  "dynamics" 


-u?- 

(u>  •  1*  is  large  enough  in  this  case, see  Fig.  6a-c) .  Note  how  this  frequency 
response  approach  is  indifferent  to  system  order. 


Step  2.  Bounds  B(u>)  on  G(jw)  Given  the  set  g  »  (P)  ,  the  problem  is 
to  find  F(s)  ,  G(s)  in  Fig.  1  such  that  the  system  transfer  function 


T(ja<)  -  FGF/(  l+GP)  6 


for  nil  P  € 


One  may  program  the  computer 


to  find  the  (unique)  hounds  B(a>)  on  G(Jw)  ,  so  that  as  P  ranges  over 

,  .Mn|r(Jio)|  *  Ain  1  1  5  (A^fm)  -  (u>)  ]  dh  of  Fig.  4.  Alternatively, 

this  may  he  done  by  hand,  giving  useful  insight:  The  template  of  l,(|«i  )  - 
P(J«Oj)  is  that  of  P,  shifted  (in  the  Nichols  chart)  by  Arg  G( J.o  ) 
in  the  x  axis  and  by  20  1  og^  |G( j..i^)  |  in  the  v  axis.  Suppose  ^Jp(|2) 
is  given  by  ABCl)  in  Fig.  7  and  one  tries  positioning  it  at  A'B'C’O'  ,  to 
give  the  template  of  L(j2)  .  From  the  contours  of  constant  |l,/(l+l,)|  in 
the  Nichols  chart ,  it  is  seen  that  the  maximum  change  in  | T ( 1 2 ) |  is  then 
closely  (-.4*1)  -  (-5.7)  »  5.2  dh  ,  with  maximum  at  C'  ,  minimum  at  A'  . 


Suppose  the  |T(.]m)|  tolerance  permits  a  maximum  change  of  6.5  dh  at 


so  the  above  trial  is  conservative.  The  template  may  bo  shifted  lower  to 
A"B"C"D"  ^  wl\ f oh  the  |T(jw) |  tolerances  arc  precisely  satisfied.  Choose 

some  specific  P  as  nominal,  e.g.,  point  A  giving  Pq(|2)  *  -11  db  /- 1 53.4* 


Since  l,_(|2)  (A"  in  Fig.  7)  ■  -  4.2  db  /  -60°  ,  the  corresponding  bound  on 


G(J 2)  -  !.0(|?)/P0(.12)  Is  8.8  db  /93.4e  ,  l.e.,  if  Arg  G()2)  is  03.4*  . 

It  Is  necessary  for  | G( j 2 ) |  *  8.8  db  ,  in  order  that  A | T ( J 2 ) |  t  b.5  db  duo 
to  the  uncertainty  In  P(J2)  .  This  manipulation  of  p(2)  is  repeated 
along  a  new  vertical  line,  giving  another  point  on  the  boundary  B(2)  of 
permissible  G(J2)  . 
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Fig.  8  shows  the  bounds  so  obtained  on  G(jm)  ,  and  the  G(ju>)  chosen 
to  satisfy  these  bounds.  Let  e^  be  the  excess  of  poles  over  zeros  assigned 

to  G(s)  ,  so  that  as  s  -*  »  ,  G  -*  k  /s  .  It  is  reasonable  to  define 

G 

the  optimum  G  as  that  which  satisfies  its  bounds  with  minimum  k(,  .  It  has 

boon  shown ^  that  G  lies  on  B(u<)  at  all  w  and  that  G  exists  and 
opt  opt 

is  unique.  The  design  of  a  practical  G(s)  to  satisfy  the  bounds  is  some- 

3 

what  of  an  art  .  For  a  given  skill  in  the  art,  the  greater  the  number  of 
poles  and  zeros  of  G  ,  the  closer  one  can  got  to  the  optimum,  so  there  is 
trade-off  between  complexity  and  bandwidth.  Hero,  wo  chose  simply  by  cut 
and  try  G(s)  *  (1  +  .2s/s(l  +  .033s  )(1  +  .002s  +10  *s  )  with  very  modest 
bandwidth,  Figs.  4,8.  A  much  simpler  G(s)  could  have  been  chosen  with 
larger  bandwidth.  The  designer  must  make  his  own  trade-off.  Ref.  3  offers 
some  advice  on  the  shaping  of  a  function  to  satisfy  a  set  (R(<»)}  in  the 
Nichols  chart. 

Step  3.  Design  of _ F(s)  ■  G(s)  only  guarantees  that  a|T(  |<o)|  <  A^(<i>)-A^(c.') 

of  Fig.  4,  o.g.,  at  hi  «  10  i  the  actual  change  in  |L( j  10) / (1+1. ( j  10)>(  is 
from  -7  db  to  4  db  ,  while  from  Fig.  4,  the  permitted  change  in 
1 T( J 10) |  =  |FL/(1+L) j  is  from  -15  to  2.8  db  .  Hence,  any  value  of 
| F( J 10) |  €  [-8,  -1.2  db]  is  acceptable.  In  this  way,  upper  and  lower  bounds 
on  In  jio)  |  are  obtained  and  F(s)  is  chosen  to  satisfy  them,  which  is  also 
somewhat  of  an  art.  In  this  example,  a  satisfactory  F(s^  -  1 1+ .  33sH  1  +  .05sl 
(1+. 25s)(l+. 2s) (1+. 0125s)2, 


see  Fig.  4. 
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Deslgn  Results.  The  nonlinear  system  was  simulated  and  its  response  found 

for  several  hundred  command  inputs  and  gust  disturbances,  some  of  which  are 

★ 

shown  in  Figs.  3,  9-12.  Some  typical  responses  to  c  step  commands  are  shown 
in  Fig.  3  for  various  ,  p  and  step  (k)  values.  The  assigned  bounds  on 

A  9 

c  are  also  shown.  The  transient  response  of  a(t)  ,  0(t)  etc.  depend,  of 
course,  on  the  values  of  k  ,  p  ,  .  Two  sets  of  these  are  shown  in  Figs. 

9a,  b,  with  Fig.  9a  depicting  very  large  a(t)  excursion,  for  which  the 
Cij  (a)  in  Fig.  2  are  in  strongly  nonlinear  ranges.  These  are  the  inputs 
for  which  the  system  was  designed,  and  for  which  guarantees  can  he  made.  In 
a  very  few  cases  there  was  very  slight  excursion  out  of  the  bounds.  It  is 
possible  to  Include  in  the  design  other  inputs  and  gust  disturbances,  with 
specified  response  tolerances  -  and  then  guarantees  can  be  made  for  these  as 
well.  The  response  to  other  inputs  is  nevertheless  found  here  to  be  also  quite 
satisfactory.  This  is  typical  of  the  design,  i.e.,  the  system  is  not  very 
sharply  tuned  to  the  class  of  inputs  used  in  the  design  execution.  There  is 
reasonable  response  continuity  to  other  inputs. 


Some  responses  to  very  large  c  step  commands  causing  hard  6  saturation 
are  shown  in  Fig.  10.  Response  to  Gust  Disturbances.  The  gust  input  was 


modelled 


by  replacing  a  in  (3^  by  a  *=  tan  —  +  a  .  .  Two  kinds  of 

VQ  gust 


a  were  used.  In  one  a  is  a  half-sine  wave  of  amplitude  ?0/V„  radians 

gust  gust  0 

and  half-period  THALF  6  [.2,2] sec  .  Some  results  are  shown  inf  Figs.  11a, b. 

In  both,  the  gust  begins  precisely  at  the  instant  of  application  of  simultaneous 

* 

c  step  commands.  The  second  kind,  as  in  Fig.  lib,  is  stochastic  guassian  with 

2  2 

power  spectrum  k/(l+w  )V„  and  (a  1  =  6/V„  radians.  Examples  of 

r  0  gust  rms  0 

k 

responses  to  a  single  square  wave  c  command  with  equal  positive  and  negative 
values  k  and  total  duration  2  THALF  are  shown  in  Fig.  12. 


DISCUSSION  AND  CONCLUSIONS 


A  ««cond-order  model  was  used  lor  T(s)  with  ex. css  of  poles  ovci  zeros 

•  •  1  ,  In  Step  1.  This  appears  to  ho  Incompatible  with  0,  F  fifth  order 

and  excesses  e„  •»  A  ,  o  -  2  ,  and  P  (no  nnalvtle.il  expression  but  linear  I  zed 
u  r 

(l  A)  give  second-order  1’)  with  Cj,  -  0  .  Strict  design  execution  appears 
to  require  a  T(s)  of  complexity  compatible  with  T  -  V  Gl*/(1  4  Cl’).  Note, 
however,  that  In  Step  1,  the  model  of  T(s)  Is  used  only  to  generate  a  set 
which  covers  the  range  of  apriorl  specified  acceptable  outputs  Any  T(s) 

model  which  achieves  this  Is  clearly  satisfactory,  and  the  simpler  the  hot  let  . 
The  designer  can  later  choose  the  complexity  of  C.(s)  ,  K(s)  ,  witli  no  regard 
for  that  of  the  T(s)  model  used  In  Step  i. 

The  class  of  applicable  non  I  Inear i t les  has  been  defined  ImplJcitly  In 
Ref.  2,  but  one  very  large  class  can  be  dc lined  explicitly.  Let  Dj  v(t)  " 

D2  x(t)  with  ,  D^  operators  which  may  be  nonlinear,  uncertain  and 
time-varying,  e.g. 

n  u/'w...5  .  .  (  A  +  lit  e-'*1  \  .?  I  in  ? 

D.  y  «  M(y)(y)  sgn  y  4  !  -  j  y  |y|  4  II  v 

'  F.  4  Ft  ' 

E  €  (1,5)  ,  F  €  ( . 5 , A 1  ,  A  €  [-1,6]  ,  a  €  (.5,1.5) 

B  €  (-3,2)  ,  H  €  ( -A . 1 )  ,  n€  (.5,2)  ,  M€  (1,5)  . 

The  range  of  M  must  be  of  the  same  sign.  All  y  €  ^  and  y  must  be 
bounded  for  all  t  In  ( 0 ,«* )  and  y  must  exist.  Hence  y  must  be  twice 
differentiable  except,  at  most,  at  a  countable  number  of  points.  Then 
Dj  y  »  |(t)  Is  known  and  there  must  exist  a  unique  solution  for  x  in 
Dj  X  •  (Kt)  .  The  solution  must  he  bounded  for  all  t  E  (0,“’)  .  Thus, 

D2  *  "  <|»(t)  must  be  "bounded- Input ,  bounded-output"  stable.  However, 

z(t)  «  v(t)  may  be  "unstable"  In  that  a  bounded  v(t)  in  allowed  to  resell 
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in  unbounded  z(t)  .  It  is  only  necessary  that  bounded  z(t)  gives  bounded 
v(t)  . 


It  is  possible  that  a  simple  linearization  might  do  just  as  well  in 
practice,  but  this  is  a  matter  of  chance,  whereas  this  design  technique  is 
guaranteed  to  work  if  the  constraints  are  satisfied.  It  can  be  argued  that 
the  constraint  of  uniqueness  may  bo  waived  in  Sec.  2,  if  the  set  of  inputs 
(x(t)}  which  can  give  the  single  output  y(t)  ,  is  "compact",  l.e.  if  that 
set  can  be  covered  as  accurately  as  desired  by  a  finite  number  of  elements. 
A  single  w  ,  y  pair  then  generates  a  set  of  Iti  P  instead  of  only  one. 
However,  this  has  not  been  rigorously  proven  as  yet. 


The  constraint  on  7h  that  is  mp,  is  required  because  only  then 

can  one  guarantee  that  any  specifications  no  matter  how  narrow  (A2-Ai 
arbitrarily  small  but  nonzero  in  Fig.  4), may  be  satisfied  for  arbitrarily 
large  hut  bounded  parameter  uncertainty  (but  some  parameters  must  not  change 
sign^) .  No  such  gurantcc  can  be  made  for  nmp  but  the  problem  is  still 

solvable  if  the  specifications  are  not  too  narrow  and  j^is  not  too  large®*^ 
a  set. 


Finally,  a  recent  extension^is  to  linear  and  nonlinear  multiple  input- 

2 

output  systems,  where  the  problem  is  transformed  into  the  design  of  n 
single-loop  ltl  systems  like  Fig.  1. 
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FIGURE  TITLES 


System  structure. 


is  replaced  by 


(p)  . 


Aerodynamic  coefficients  C^(a)  . 

Normalized  bounds  on  acceptable  c*(t)  .  Representative 
simulation  results. 


Bounds  on  | T ( ju>)  |  .  Bode  plots  of  designed  |F(ju>)|  ,  |G(j<,»)| 
Loci  of  typical  Iti  P(jw)  .  Cases  e,f  are  open-loop  unstable. 
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c.  Templates  c/p(w)  of  {P(ju>))  at  u<  =  .04,  1.9,  1?  . 


Hand  derivation  of  bounds  on  L(jw)  in  Nichols  chart. 


Bounds  B(w)  on  G(jw)  ,  and  G(jut)  chosen. 


i 


Representative  responses  of  a(t),  e,  6,  c*. 

Responses  for  Input  causing  hard  <5  saturation. 

Responses  to  gusts  and  simultaneous  c*  step  command 
(a)  half-sine  gust,  (b)  random  gausslan. 

Responses  to  square  wave  c*  command. 
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FUNCTIONS  IN  DEGREES 


APPENDIX  2 


Derivation  of  the  transfer  function  of  a  system. 


from  an  input -out put  pair 


Let  x(t)  ,  y(t)  be  the  input  and  output  of  a  system  for  t  6  [ 0 , T )  . 

Find  a  transfer  function  T(s)  -  Y(s)/X(s)  ,  where  Y(s)  »  Ly(t)  ,  etc. 

Actually,  we  seek  a  rational  approximation  of  T(s)  ,  P(s)/Q(s)  ,  where  P 

m  . 

and  Q  are  polynomials  of  degree  m  ,  n  respectively.  Let  P(s)  -  £  Pc*  • 


Q(S)  '  1  q"=  ’  * 
k-0 


SM  y(s)  -  X (s ) 


Y  (s)  P(s)  Q(s)Y(s)  -  P(s)X(s)  Sn  Sr 

xfSTo&  xTTjQjs)  QTiTx(sy 


If  VTTV  is  a  rational  function  then  —-S—  y(s)  -  ^ - S •  x(s)  *  0  is  an 

'  sn  sn 

exact  representation.  We  seek  the  coefficients  pk  ,  q^  which  minimize 


Q(  jw)  Y(  jto) 
(jw)" 


( j  w)  0 


n  qkY(jw)  m  p  X(jw)i2 

'  1  A  77^ '  A  7^*  d“ 

*®  k-0  (Jw)  k-0  (jw)  1 

00  7 

Parseval's  formula  gives  J  -  j  |f(t)|  dt  .  where 


n  /  i  \  fD 

f(t)  -  l  qky'(n‘k)(t)  -  l  p 
k-o  k-0 


x-(n‘k)(t)  . 


/_ i \  t  t 

y '  (t)  ,  the  l-times  integral  of  y(t)  ,  J  ...  /  y(t)d(.‘  . 

0  0 

J-  /  |f(t)|2dt  +  /  |f  (t)  |2dt  ♦  /"|f(t)|2dt  ;  /  |f  (t)  |2dt  -  o  . 

-o’  0  T 

since  y(t)  -  0  ,  x(t)  -  0  for  (.0  ;  /  jf(t)  |2dt  -  0  ,  since  x(t)  , 

T 

y  ( t )  for  t  }  T  can  be  defined  correspondingly  for  any  p^  ,  .  Therefore 


H 
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J  »  /T|  f  (t)  | 2  dt  l1*) 

0 

and  the  minimization  of  J  determines  uniquely. 

Let  x(t)  ,  y (t)  be  given  at  the  points  t = 0 ,h ,2h , . . . Ih  -  T  .  We 
construct  the  integral  y^  ^(t)  by  the  following  formulae: 

y(l) (0)  =  0  ,  y(_1> (ih)  =  y(_1) (h(i-1))  +  j-  |y( ih)  + y(i-l)h) ] 

The  multiple  integrals  y^  ^  (t)  ,  U  ks  n  and  x^  U^(t)  ,  n-m^u^n  are 
computed  in  the  same  manner.  We  replace  the  integral  (k)  by  the  corresponding 
sum: 

,  .  v  T  V  (n-k),.^  v  (n-k) / ..  \  ]2  , 

J  =  h  I  I  <l|,y  (ih)  -  I  p.x  (ih)  ,  q  -1 

i=0  L  k=0  k=0  J 

and  determine  the  coefficients  p^  ,  q^  as  coefficients  of  linear  regression 
(by  the  least  square  method). 


Remarks 

(1)  The  degrees  of  the  polynomials  P(s)  and  d(s)  should  be  chosen 

from  the  apriori  information  about  the  system.  If  apriori  information  is 

lacking,  different  degree  choices  of  P  ,  Q  may  be  tried  and  the  best 
selected.  Another  kind  of  apriori  information  which  could  be  used  in  this 
method  is,  for  example,  asymptotic  behavior  of  the  functions.  The  difference 
in  degree  of  P  ,  Q  can  be  judged  by  the  difference  in  smoothness,  at  t*0  , 
of  x(t)  and  y(t)  ,  because  of  the  initial  value  theorem  in  Laplace 
transforms. 

(2)  If  the  function  P(s)/Q(s)  is  stable,  then  the  approximation  can 

be  checked  by  computation  of  y(t)  =  [ x ( t )  I  and  comparison  of  y ( t ) 

with  y(t)  . 
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